This notebook demonstrates how to compute theoretical bounds on RV-precision given a synthetic spectrum and an LSF. Korg has a few functions that make this easier, based on Bouchy et al. 2001, A&A, 374, 733.

In addition to Korg, we'll use PythonPlot in this example. If you don't have it installed, you can run using Pkg; Pkg.add("PythonPlot") to install it.

    using Korg, PythonPlot

First, specify the wavelengths and line spread function (LSF). We'll use the APOGEE instument's, wavelength coverage and resolution for this example. Reduced APOGEE spectra are resampled onto wavelengths which are uniform in log-wavelength.

delLog = 6e-6;
apowls = 10 .^ range((start = 4.179 - 125 * delLog); step=delLog, length=8575 + 125)
8700-element Vector{Float64}:
 15074.74588598709
 15074.954152736687
 15075.162422363615
 15075.370694867916
 15075.578970249626
 15075.787248508786
 15075.99552964544
 15076.20381365962
 15076.412100551372
 15076.620390320733
     ⋮
 16997.928574606198
 16998.163411290072
 16998.398251218357
 16998.633094391098
 16998.86794080834
 16999.10279047013
 16999.33764337654
 16999.572499527556
 16999.807358923248

We compute the LSF from the APOGEE resolution ($R \approx 22,500$) as a sparse matrix. (See the LSF explanation here for more details.)

LSF = Korg.compute_LSF_matrix((15_000, 17_000), apowls, 22_500)
8700×200001 LinearAlgebra.Adjoint{Float64, SparseArrays.SparseMatrixCSC{Float64, Int64}} with 2103823 stored entries:
⎡⠀⠈⠉⠉⠓⠒⠒⠒⠒⠶⠤⠤⠤⠤⢤⣀⣀⣀⣀⣀⡀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⎤
⎣⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠉⠉⠉⠉⠉⠙⠒⠒⠒⠒⠒⠦⠤⠤⠤⠤⢤⣀⣀⣀⎦

Now, let's synthesize a spectrum, and simultaneously apply the LSF and resample to the APOGEE wavelength grid.

apolines = Korg.get_APOGEE_DR17_linelist(; include_water=true)
wls, flux, _ = Korg.synth(; linelist=apolines,
                          wavelengths=(15_000, 17_000),
                          Teff=5777,
                          logg=4.44,
                          M_H=-1.1,)

SNR = 50.0
obs_err = (LSF * flux) ./ SNR;

Exact radial velocity precision

This function provides the exact RV precision bound (in the small redshift limit) in m/s. Note that we are using the unconvolved flux here (flux), not the values with the LSF applied (LSF * flux).

noise_prec = Korg.RV_prec_from_noise(flux, (15_000, 17_000), apowls, LSF, obs_err)
146.81420525789295

$Q$-factor calculation

We can also calculate the $Q$-factor, which provides an approximately $S/N$-independent notion of the radial velocity "information" in the spectrum.

Q = Korg.Qfactor(flux, (15_000, 17_000), apowls, LSF)
413.70095102641835

We can calculate the appriximate RV precision bound from $Q$. Technically, the precision depends on the root-mean-squared per-pixel $S/N$ value, but for non-line-blanket spectra, this is approximately equal to the $S/N$. (Verifying this fact is left as an excercise to the reader.)

SNR_RMS = SNR # approximately
Npixels = length(apowls)
Q_prec = Korg.RV_prec_from_Q(Q, SNR, Npixels)
155.38337597445278

This is a few m/s off of the exact calculation, but it's pretty close!

$Q$-factor as a function of metallicity

Let's examine how the $Q$-factor changes as we change the metallicity for a star with solar $T_\mathrm{eff}$ and $\log g$. As we add more and deeper lines, there is more for an RV measurement to grab onto.

M_Hs = -5:1.0:0 #metallicity ([M/H]) values
Qs = map(M_Hs) do M_H
    _, flux, _ = Korg.synth(; linelist=apolines, wavelengths=(15_000, 17_000), Teff=5777, logg=4.44,
                            M_H=M_H,)
    Korg.Qfactor(flux, (15_000, 17_000), apowls, LSF)
end

scatter(M_Hs, Qs)
xlabel("[M/H]")
ylabel(L"$Q$ [m/s]")
Example block output

We can use the $Q$-factors to calculate the RV precision across $S/N$. Note that these precision estimates are not the whole story at high $S/N$ because effects like granulation become important.

SNRs = 1:1:1000
for (M_H, Q) in zip(M_Hs, Qs)
    plot(SNRs, Korg.RV_prec_from_Q.(Q, SNRs, Npixels); label="[M/H] = $M_H")
end
legend()
ylabel("RV precision [m/s]")
xlabel("S/N")
yscale("log")
xscale("log")
Example block output